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The static and dynamic properties of many-body quantum systems are often well described by collective 
excitations, known as quasiparticles. Engineered quantum systems offer the opportunity to study such emergent 
phenomena in a precisely controlled and otherwise inaccessible way. We present a spectroscopic technique to 
study artificial quantum matter and use it for characterizing quasiparticles in a many-body system of trapped 
atomic ions. Our approach is to excite combinations of the system’s fundamental quasiparticle eigenmodes, 
given by delocalised spin waves. By observing the dynamical response to superpositions of such eigenmodes, we 
extract the system dispersion relation, magnetic order, and even detect signatures of quasiparticle interactions. 

Our technique is not limited to trapped ions, and it is suitable for verifying quantum simulators by tuning them 
into regimes where the collective excitations have a simple form. 


PACS numbers: 03.65.Ud, 75.10.Pq, 37.10.Ty, 37.10.Vz 

The level of experimental control over engineered many- 
body quantum systems has been rapidly improving in recent 
years. Atoms in optical lattices [1] and ions in electrodynamic 
traps [2], for example, now offer unprecedented access to the 
states, dynamics, and observables of interacting quantum sys¬ 
tems. As these engineered systems become larger and more 
complex, conventional tomographic methods [3] for charac¬ 
terising their states and processes soon become either ineffi¬ 
cient or of limited applicability [4]. New methods are required 
to verify what has been built in the laboratory and to study the 
properties of these fascinating systems. 

In their low-energy sector, interacting many-body systems 
are often well described in terms of collective excitations, with 
effective mass, dispersion relation, and scattering properties. 
Such emergent excitations can be understood as weakly inter¬ 
acting quasiparticles, which are responsible for a range of dy¬ 
namical properties of the underlying system [5], Recently, ex¬ 
periments in systems of atoms and ions have demonstrated the 
central role that quasiparticles play in the transport of informa¬ 
tion and entanglement [6-8]. The ability to precisely measure 
the properties of a system’s quasiparticles, therefore, becomes 
a valuable tool for studying many-body quantum systems in 
the laboratory. 

Spectroscopy is an established approach to studying nat¬ 
ural quantum systems and their emergent phenomena, well- 
established techniques including, e.g., photoemission spec¬ 
troscopy in electronic systems [9] and neutron scattering in 
magnetic materials [10]. Currently, there is important theoret¬ 
ical [11-14] and experimental [15-23] progress on developing 
spectroscopic techniques for engineered quantum systems. In 
this Letter, we present a spectroscopic technique for charac¬ 
terizing the low-lying energy spectrum of an engineered quan¬ 
tum many-body system, and apply it to study emergent quasi¬ 
particles in a system of trapped atomic ions. Our approach ex¬ 
ploits single-particle (-ion) control to excite individual quasi¬ 
particle modes. Measuring the system’s dynamical response 
to superpositions of such single excitations, allows the abso¬ 



Figure 1: (a) Energy spectrum of model Hamiltonian Eq. (1), for 

B » J and 7 spins. The degeneracy of subspaces with different exci¬ 
tation number n is lifted by the spin-spin interactions J. The system 
response to superpositions |0) + | k) (blue spheres) provides informa¬ 
tion about the absolute energy E k of spin-wave state \k). The response 
to superpositions \k) + | k') (red spheres) enables high-precision mea¬ 
surements of the system dispersion relation. Quantum dynamics in 
higher excitation manifolds (n > 2) involve both spin-wave disper¬ 
sion and scattering, (b) Spin-wave amplitudes A k .. 


lute quasiparticle energies and system dispersion relation to be 
determined. Furthermore, by probing the response to super¬ 
positions of multi-quasiparticle states we are able to resolve 
shifts in the energy spectrum due to quasiparticle interactions. 

Our system is a ID chain of 40 Ca + ions in a linear ion 
trap. In each ion j = 1 ... N, two long-lived electronic states 
represent a quantum spin-1/2 particle with associated Pauli - 
matrices <x^ (J3 = x,y,z). Under the influence of laser-induced 
forces, the system is ideally described by a one-dimensional 
transverse Ising model [24, 25], with Hamiltonian 

H - // X Jipr'jiT) + TiB X o^j , d) 

i<j J 

where cr * = (crj+Zcrp are spin raising and lowering operators. 
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The spin-spin coupling matrix J, ; ~ J/\i - j\ a approximately 
depends as a power-law on spin separation |; - j |, where 0 < 
a < 3 [8, 26]. All experiments presented here are carried out 
in a system of N = 7 spins where a ~ 1.1, corresponding to 
interactions of intermediate range [27]. 

We work with a strong transverse field B » J, in which 
case the energy spectrum of H splits into N + 1 subspaces, 
each containing eigenstates with the same integer ‘excitation 
number’ n = ^((cH) + l)/2, i.e., total number of up-pointing 
spins, t (see Fig. 1(a) and [25]). In this limit, H conserves n. 
Therefore the subspaces are uncoupled and there is no popula¬ 
tion transfer between them during time dynamics. The eigen¬ 
states in the single-excitation subspace (n = 1) are delocalised 
standing waves of spin excitation | k) = A k IT;)®li)' 8(iV_1) , 
illustrated in Fig. 1(b). These states | k) are the fundamental 
quasiparticle modes of the system, with eigenenergies E k that 
define the dispersion relationship. 

At the basis of our spectroscopic method is the preparation 
of superpositions of eigenstates. Specifically, beginning with 
the ground state |J,) 0jv = |0), we use a laser, focused to a sin¬ 
gle ion and rapidly switchable between different ions, to ro¬ 
tate each spin j into the state | Of) = cos (Oj'j |J, j) + sin \0jj IT j). 
As we will show, preparing different superpositions allows us 
to spectroscopically measure energy gaps from the beat notes 
observed between eigenstates. 

Spectroscopy of ground-state gaps .— First, we measure the 
energies E k of the single-excitation eigenstates | k) relative to 
the absolute ground state. Although the | k) modes are delo¬ 
calised, entangled states, we can prepare product states that 
are good approximations by setting 8j = tan -1 (yAj), yielding 
the /V-spin input state 

N 

\<ffk) « 0 [li j) + jA k j I T;>] = |0> + y | k) + O (r) . (2) 

j= i 

Preparing this state requires no detailed a priori knowledge 
of the Hamiltonian. As we will see, a strong signal is 
obtained when modelling the exact eigenmodes by generic 

standing waves with amplitudes A k = y sin For 

open boundary conditions and nearest-neighbour interactions 
(a = oo), these form the exact eigenstates. In case of finite- 
range interactions (a < oo), the amplitudes A k are modified, 
but the number of wave-function nodes remains an accurate 
way to characterize the eigenstates. 

After preparing the system in \f k ), we let it evolve to 
| f{t)) = exp (-iHt/ti) |t (fk), and measure its single-spin magne¬ 
tizations in the x-y plane of the Bloch sphere. Specifically, we 
measure in a rotating frame where x = cos(2 Bt)x - sin(2/ir)y 
and y = sin(2 Bt)x + cos(2 Bt)y [25], The corresponding 
expectation values are {(T x j y \t)) = {if/(t)\ cr x j y) \f(t)) ~ 

2y Re[exp(-iskt/fi) (0\ |k)] + C?(y 3 ), where h = h/ln is 

the reduced Planck constant. These expectation values oscil¬ 
late at a frequency determined by the energy shift Sj, = E k -2B 
that result from the spin-spin couplings J,j [Fig. 1(a)]. The 
shift Ek can then be extracted via a Fourier transform. To 
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Figure 2: Energy measurements of quasiparticle states, (a-c) An 
absolute energy probed with initial state |0) + \k = 7) and y = 0.4. 
(a) Spin-resolved x-magnetization dynamics, (b) Summed signals 
Wjjjft) along x and y. (c) The Fourier transform of MAt) + iMy(t) 
shows a peak at -249±42 Hz, revealing the downward energy shift of 
the |k = 7) state caused by the spin-spin interactions, (d-f) A relative 
energy probed with initial state \k = 1) + \k = 7) and y = 0.7. (d) z- 
magnetization dynamics, obtained by post-selecting outcomes with 
exactly one spin-excitation, (e) Summed signal Mft), see text, (f) 
The Fourier transform of M,(t) shows a peak at 383± 13 Hz, revealing 
the energy splitting of the two quasiparticle states. In (a-c) y = 0.4, 
corresponding to a probability of having zero (one) excitation of p 0 = 
0.74 (pi = 0.24). For (d-f) y = 0.7, yielding po = 0.43 and p\ = 
0.42. 


maximize the signal, we transform Mft) + iMft), where 
Mu (y)(t) = X/ A k j{tr x - ik> (t)), and mirror the result to negative 
times to minimize transform artifacts. For these measure¬ 
ments, we choose y such that the probability of creating more 
than one quasi-particle is only p 2 = 0.02, thus avoiding spu¬ 
rious frequency components in the Fourier transform due to 
energy gaps between states with one and two quasiparticles. 

Results for \fk=i) are presented in Fig. 2(a-c). The Fourier 
transform reveals a distinct peak at Ek=i/h = -249 ± 42 Hz, 
which reveals the interaction-induced dispersion s ,t of mode 
| k) around the band center 2 B. From the observation that 
Ek = 7 < 0 and the anti-aligned nearest-neighbour amplitudes 
of | k = 7), we can infer that J t j > 0, i.e., the spin-spin cou¬ 
plings are antiferromagnetic. Analogous results are found for 
other modes. We observe, e.g., a positive energy shift for E k -.i, 
in agreement with J, j > 0 and the requirement e k = 0. It 
is worth noting that the measured quantity M^(t) is equal to 
the quasiparticle Green’s function Gift), which contains all es¬ 
sential information about the quasiparticles, such as energies, 
scattering rates, and life times [11, 28], 

Spectroscopy of the dispersion relation. — The measure¬ 
ment precision of the absolute energies E k is limited by de¬ 
phasing due to laser-frequency and ambient magnetic-field 
noise. As subspaces with fixed excitation number n are 
decoherence-free with respect to these noise sources [8], the 
relative energies A E k ,k' = Ek - can be measured with 
much higher resolution by preparing superpositions of the 
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Figure 3: Quasiparticle dispersion relation, obtained from Fourier 
transforming the measured dynamics of the initial states \k = 1) +1 k') 
with k' = 2... 7 (blue traces). The blue circles mark the positions of 
peak maxima which closely match the theoretical dispersion relation 
(red line). 


type | k) + | k'). 

Here, we initialize the system in a product state similar to 
Eq. (2), and choose the rotation angle Oj = tan -1 | y (A k + Aj)] 
to create the state 

Wkk') = |0> +y{\k) + \k')) + 0(y 2 ). (3) 

We then measure the time-evolved state in the z-basis. Since 
H commutes with such measurements, this allows post¬ 
selection of the | k) + | k') contribution to the initial state 
by retaining only measurements projected on the single¬ 
excitation subspace. The individual z-magnetizations are 
added up to obtain the signal M\(t) = C/(o~(?)) = 

cos(AE c jiklcr-jlk’) with cj = sign [A k A k .^. A Fourier 

transform of M\(t) yields the energy difference AE k ^ between 
quasiparticle states | k) and | k'). 

Figure 2(d-f) shows experimental data for a superposition 
of the quasiparticle states with lowest ( k' = 7) and highest 
(k = 1) energies. Compared to Figure 2(a-c), the much longer 
coherence times within the n = 1 subspace result in a sig¬ 
nificantly higher spectral resolution. By repeating this exper¬ 
iment for superpositions of different spin waves, we recon¬ 
struct the complete single-quasiparticle dispersion relation, 
measured relative to k = 1. In Fig. 3, the results are com¬ 
pared with theoretical predictions, obtained by diagonaliza- 
tion of the model’s coupling matrix Jjj, with J as only free pa¬ 
rameter. Instead of a cosine dispersion typical for short-range 
physics, the data reproduces very well the large group velocity 
(da)/dk) at k = 1, as expected for long-range interactions with 
a < 2 [27, 29]. The fact that each of the N frequency spec¬ 
tra in Fig. 3 is dominated by a single peak demonstrates that 
the heuristically chosen sine-wave coefficients A k . are, indeed, 
close to the ones of the exact eigenstates. 

Spectroscopy of interacting quasiparticles .— 

Until now, we probed single quasiparticle dynamics. In 
the spectroscopy scheme based on initial states like Eq. (3), 


quasiparticle excitation numbers above n — 1 are naturally 
realized in higher orders of y. We can therefore go further 
and perform spectroscopy in higher-excitation sectors by post- 
selecting time-evolved states with n spin excitations. This 
provides an opportunity to characterize interactions between 
pairs of quasiparticle excitations. 

We now prepare product states \f k k f, Eq. (3), with 

’ . ( 2 ) 

larger y, so that the next-order contribution W kk ,) = 

7 2 [2 j + Aj) cr|] 2 = y 2 (| kk) + 2 | kk') + | k'k')) becomes 
important. After evolving | ifkk') under H as before, we per¬ 
form a projective measurement of <rj(f). Afterwards, we post¬ 
select measurement outcomes with two excitations, compute 
the observables Vf = j((l + crj)(l + <xj.)), and Fourier- 

transform the sum = T,i±j (t) with suitable 

weights c" or cf ; [25]. In parallel, we post-select outcomes 
with a single excitation and thereby measure M\(t) to obtain 
the single-particle energy splitting \AE kk '\/h as a reference 
[Fig. 4(c)]. 

Without interactions between quasiparticles, higher ex¬ 
cited states would be simple combinations, determined by the 
quasiparticle statistics, of the single-excitation states. The 
measured frequency spectrum for | f k , =7 ), however, reveals 
that the quasiparticles do interact: the peak in Fig. 4(a) at 
665 + 14 Hz is significantly lower than 2\AEk,k'\lh = 765 + 
32 Hz, which would be expected for non-interacting bosons 
(NIB) from double population of the states k and k'. Further¬ 
more, the aforementioned peak is higher than the maximum 
allowed for non-interacting fermions (NIF), which is, due to 
Pauli exclusion, \E k =i + E k= 6 - E k= i - E k =i\lh = 585 ± 21 Hz. 
These measurements show that the quasiparticles in our sys¬ 
tem cannot be described by a simple free particle model, 
which is a consequence of the long-range nature of our in¬ 
teractions. A NIF model would be accurate only for strictly 
nearest-neighbour interactions [30], 

For a qualitative understanding of the two-particle spectra, 
it is illustrative to map the spin operators to hard-core bosonic 
creation and annihilation operators [5], crt —> b\ (crj —> bj ), 

with quasiparticle modes | k) = b\ |0) = |0). Note that 

the two-particle states | kk’) = b\b^ k , |0), k ( = 1... N, are not 
eigenstates of the n — 2 subspace; there would be N(N + l)/2 
different states | kk'), whereas the n — 2 subspace has only 
N(N - l)/2 dimensions because of the hard-core constraint 
imposed by having a spin-1/2 system. Nevertheless, for N » 
2, we expect | kk') to have a significant overlap with a two- 
quasiparticle eigenstate. 

Without the hard-core constraint, \wf,) would populate 
three exact eigenstates, whose energies would be the sums of 
the individual quasiparticle energies E k and E k f [see Fig. 4(d)], 
In that case, we would measure two distinct spectral peaks at 
frequencies \AE kk ’\/h and 2\AE k ^\/h. However, the hard-core 
interactions induce an interaction shift away from this value. 
For low quasiparticle densities, n N, qualitative insights 
into interaction effects can be obtained by treating the hard¬ 
core interaction as a small perturbation. First-order perturba- 
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Figure 4: Spectroscopy of quasiparticle interactions. (a,b) Fourier 
spectra of response to the two-quasiparticle state |</f®, =| 7 ) with y = 
1.4 (red curves), obtained from summed signals M 2a (t) and Ir¬ 
respectively. [25]. The highest frequency in (a) is inconsistent with 
the predictions for both non-interacting fermions (NIF) and non¬ 
interacting bosons (NIB), indicated by light, respectively dark, green 
lines. For NIB, e.g., one would expect a signal at twice the single¬ 
particle energy splitting Aplotted in panel (c). The measured 
peak positions in (a,b) compare well with those of a longer simu¬ 
lated evolution, indicated by the narrower gray curves. A pertur¬ 
bative treatment of quasiparticle interactions (blue solid half-lines) 
agrees qualitatively with the observed interaction shifts, but overes¬ 
timates the effects. Here, y = 1.4 yielding p b = 0.08, p\ = 0.33. and 
p 2 = 0.41 in order to maximize the two-quasiparticle population, 
(d) Schematic energy diagram. Left: In a non-interacting system of 
free bosons, all higher excitation eigenstates can be inferred from 
the single-particle subspace by successive application of creation op¬ 
erators. Right: For interacting hard-core bosons, the energy levels 
acquire a density- and mode-dependent interaction shift. 

tion theory predicts energy shifts of the states \kk') [25], 
and, consequently, beatnotes at v a = |A E k ,k' + Vk,k ~ V b ,k'\/h, 
v b = |A E kJc r - V k ',k' + VkjA/h, and v c = v a + v b . In Fig. 4(a-b), 
the frequency spectrum for different summations of the data is 
shown for W kk , = \ 7 ), where, indeed, we observe three distinct 
peaks near v a , v b and v c . 

While the predictions of perturbation theory prove qualita¬ 
tively correct, in terms of the existence and direction of inter¬ 
action shifts, the effects are overestimated. In addition to the 
predicted three-peak structure, and their sum and difference 
frequencies, another large peak is observed at a very low fre¬ 
quency [Fig. 4(a)], This peak can be explained by studying 
the exact eigenstates, where we find that the initial state has 
substantial overlap with a fourth eigenstate not expected from 
first-order perturbation theory. The reason for these discrep¬ 
ancies is that for the small chain considered the quasiparticle 
density, n/N = 2/7, is not sufficiently low to quantitatively 
treat the hard-core interaction as a small perturbation. With 
increasing system size n/N —» 0, the quality of perturbation 
theory will improve and the states \kk') should provide closer 
approximations to the true eigenstates of the two-excitation 
subspace of H. 

One can use this scheme to study interaction effects be¬ 
tween different combinations of quasiparticle states. The re¬ 
sulting spectra can be much more complex than the one pre- 


(2) 

sented in Fig. 4. For example, w k l = \ 4 ) has considerable over¬ 
lap with many of the exact two-excitation eigenstates, which 
can be interpreted as a strong scattering of quasiparticles into 
other modes [25], This scattering reduces the quasiparticle 
lifetime, and induces broader spectral features. 

Conclusion .— We realized a Ramsey-like technique for 
measuring the low-lying energy spectrum of an artificial quan¬ 
tum many-body system of trapped ions. Our technique en¬ 
abled a complete characterization of the system’s quasiparti¬ 
cles, including the dispersion relationship and scattering prop¬ 
erties. While the former determines how information is trans¬ 
ported in our system, the latter shows that our system can¬ 
not be described by a simple model of non-interacting parti¬ 
cles. The technique works by tuning the system into regimes 
where collective excitations have a simple form, allowing for 
many-body couplings to be precisely studied. Future work 
could study the regimes where B and J are of similar size. 
Here, the quasiparticle peak in the spin-wave Green’s func¬ 
tion would broaden and eventually disappear at the quantum 
critical point. Furthermore, it would be interesting to investi¬ 
gate whether multi-dimensional spectroscopy techniques de¬ 
veloped forNMR or femtosecond spectroscopy [32, 33] could 
be combined [34] with the single-site addressability available 
in many artificial systems to measure finite quasiparticle life¬ 
times due to their interactions or incoherent coupling to the 
environment. 
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SUPPLEMENTAL MATERIAL 

Derivation of a transverse Ising model 

This section provides a short derivation of how to real¬ 
ize a transverse field Ising model in a linear string of ions. 
The mechanism is based on transversally illuminating the ion 
string by a bichromatic laser field that off-resonantly couples 
to the lower and upper vibrational sidebands of the transition 
between the two electronic states encoding the pseudo-spin. 

We start by considering a monochromatic running-wave 
laser field coupling the internal electronic states of an ion to 
the 2N transverse motional modes of the N ion string. We 
work in an interaction picture with respect to the Hamiltonian 
Ho - ^ 2,cr: + 2m tioj m a,„a m [1], where Tico o is the energy 
splitting between the two pseudospin states and Tiu> m the en¬ 
ergy of phonon mode m with annihilation (creation) operator 
a m (aj„). The interaction Hamiltonian Hj for a monochromatic 
light field reads 

ti N 2N 

Hi = ^22 Q ‘j\ (T+ j e ~ <et ~ <S>L) ( ! + ir]j.m(a m e~ ,0)mt + al,e ,w '"')) 

7=1 m= 1 

+h.c. j, (SI) 

where flj is the Rabi-frequency of the laser at the i-th ion, 
d*/ the laser phase, e - a>i - a >o the detuning of the laser 
frequency from the atomic transition, and r; j m is the Lamb- 
Dicke factor, which is proportional to the displacement of the 
j th ion in the m ,h collective mode. The summations run over 
all N ions and over all 2N transverse modes. The spin-flip 
operators are given by cr* = <T x ±icTj, where cr*’ y ’ z are the usual 
Pauli matrices. Note that Eqn. (SI) is valid under the usual 
Lamb-Dicke approximation (r/j m <s 1), which is appropriate 
for our experiments. 

In order to achieve an effective spin-spin interaction, the 
red and the blue sidebands of all 2 N modes are off-resonantly 
driven with a bichromatic light-field with symmetric detun¬ 
ings e± — ±(u>m + A). Here, A denotes the detuning of the 
laser with respect to the blue sideband of the motional mode 
with the highest frequency to M . If we perform a rotating-wave 
approximation to retain only the coupling of the frequency 
component with e+ (e_) to the upper (lower) motional side¬ 
bands, the resulting Hamiltonian is given by 

h N 2N 

H hlc = T 2 XI T] i^i + (T ~i) ( a ’” e,A "‘' + a li e ~ ,Am ‘) ( S2 ) 

J=1 m= 1 

where A„, = A + (com - a>m)- In the limit of \rjj m Q.j\ <s A„„ the 
motional states are only virtually excited and serve to mediate 
an Ising interaction between the spins. This effective spin- 
spin interaction can be calculated from second-order pertur¬ 
bation theory, resulting in the effective Hamiltonian 

"eir = n Yj J U « + (T ~i ) H +<r l) = n Yj J,] °1 (S3) 

hj hj 
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with spin-spin coupling matrix 


Jij ~ 




2 N 

Z 


j,r 


4 ^ A 

m= 1 


^>M ~ Wn 


(S4) 


Further details of the derivation can be found in references 
such as [2] or [3]. 

If the detuning of the bichromatic laser beam contains a 
slight asymmetry, e± = +{u>m + A) + 6, S <K A, the raising 
and lowering operators appearing in Eqs. (S2), (S3) obtain a 
time-dependent phase factor, cr + —> cr + e~" 5 ' and cr~ —> <r e ,St , 
and thus 

Heff = J, > (°? a+ i e ~ 2 ' 6 ‘ + cr t°"] + o-J <r) + crj cr- <? 2 '*) . 

Uj 

(S5) 

Additional terms oc (T z .a m a m that arise from off-resonant ex¬ 
citation of the sideband transitions are smaller by a factor 6/A 
and can be neglected in the limit 6 <sc A in which we work 
here. 

The time-dependent Hamiltonian in Eq. (S5) can be trans¬ 
formed into a time-independent Hamiltonian by writing H e ff = 
(H e g - //,,) + //„ and changing into an interaction picture with 
respect to //„ = -h Z, f cH, generated by the unitary operator 

U s = e~ iHitin , 


K« = ul(H eS -H s )Us 

= h Y Jjj + TiB Y o~) - Rising ■ (S6) 

i,i j 

This Hamiltonian describes an Ising model in the effective 
transverse field B = - f. 

Our experiments are carried out in the regime 6 » 
max,j{|/y|}, i.e., that of large transverse field. In that case, the 
phases of the double spin-flip terms <x + cr + and <r cr in Eqn. 
(S5) are fast rotating and can be eliminated. Finally, moving 
into the same interaction picture as that used for Eq. (S6) one 
recovers the XY-Hamiltonian 


Kff = J,j ( <T ^ cr f + hx ) + 1,11 X °1 = Hxy ' ( S7 ) 

i,j j 


Spectroscopic signals for single quasiparticle states 

We perform spectroscopy of a trapped-ion system modelled 
by the XY Hamiltonian, Eq. (S7), of N coupled spins. As dis¬ 
cussed in Fig. 1 in the main text, the spectrum of Hxy splits 
into uncoupled subspaces with integer total excitation number 
(number of up-pointing spins), allowing us to study each sub¬ 
space individually or to compare the energy of states within 
one subspace against the energy of a reference state in another 
subspace. In the single-excitation subspace, one can calcu¬ 
late the eigenenergies (£)-) and eigenmodes (\k) = cr+ |0), with 


erf = Z i A k .a z ) from the model Hamiltonian simply by diag- 
onalizing the N xN matrix 7,,. 

To spectroscopically measure the energies £)., we use su¬ 
perpositions of eigenstates. These can be approximated by 
initial product states like 

N 

Wo) = §§ |cos(0 ; ) ||> ; + sin(0 ; ) |T> ; | (S8a) 

j= i 

- cjl + Y t a n(0/)cr| + Y tan(6»,■) tan(£» J )cr ! + crj 

j Uj 

+<9(tan(0) 3 )] |0> . 

Here, C = fl 7 i cos (6j) provides the overall normalization of 
the wave function. 


Absolute quasiparticle energies 


For the spectroscopy of absolute energy values Ek, we 
choose tan (Of) — yA k to generate the initial state 


I i/a> = C 


io> + y Yi Ak j a . 1 io> + ^(y 2 ) 


(S9) 


The time-evolved state is then \fi(t)) = q~ iHxy, \ if/ k ). If this state 
is measured directly in the v-basis, with respect to the frame 
of Eqn. (S7), the signal would read: 

= (tykl e iHxrt of Q~ iHxrt \(f/ k ) (S10a) 

= C 2 y(e- iEt,IR <0| of erf |0> + c.c.) + 0(y 3 ) 
= C 2 y(e~ iEtt/n A k + c.c.) + 0(y 3 ) 

This observable contains a strong oscillation at the quasiparti¬ 
cle energy Ek. The contributions from individual spins (ions) 
can be constructively summed to Zi A k (f//(t)\ cr v \4/(t)) = 
(0| crf(t)crf |0) + c.c. + 0(y 3 ), with o~ k (t) = e' Hxyt erf e~ lHxyt . 
Note, to leading order this is equivalent to the single-particle 
Greens function [4], 

However, we do not measure directly in the x-basis (or y- 
basis) of Eq. (S7), as assumed in Eq. (S10). Measurements 
of the magnetization in the x-y plane of the Bloch sphere are 
carried out by first performing n/2 pulses with a laser res¬ 
onant with the pseudo-spin transition followed by a fluores¬ 
cence measurement projecting the spin onto either | f) or | X). 
This corresponds to measurement of (cr x ) or (cr v ) in the ref¬ 
erence frame of the Hamiltonian (S2). When changing into 
the reference frame of the Hamiltonians (S6) or (S7), these 
single-spin measurements correspond to measurements of 

cr' = e^o-'e- 1 ^ = cos(*)cr l - sin(bf)o- v (Sll) 

o-f = e l ^o- y e-‘^ = sm(6t)cr x + cos(6t)cr y (S12) 

a 2 = e‘^cr z e^ = cH . (S13) 
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Note that measurements in the logical (z-basis) are unaf¬ 
fected by this frame change. The consequence of measuring 
in a rotated frame in the x or y basis is simply that we observe 
oscillations due to energy gaps e k = E k - 2 B, as presented in 
Fig.2 (a-c) in the main text. 

Relative quasiparticle energies 

To observe beatnotes between different eigenenergies 
within the single-excitation subspace, we prepare approx¬ 


imate superpositions of two eigenfunctions, tan(0,) = 

y (A k ' + A^ 2 ), yielding the initial state 


I<Ao> = C [l + y(cr+ + o-l ) + y 2 {o- + h + <r+) 2 + <9(y 3 )] |0) . 

(S14) 

When measuring the observable cr: for spin i in the time 
evolved state \ip(.t)) = e~‘ Hxrt \t//o), the observed signal is 


<m\ oi \m> = m e iHxrl ^ e - iH ^ i^ 0 > = 


c 2 


<0| cr] |0) 


+y 2 <0| (o- k e iEt > ,,h + o-- k Q iE ^ tlh )(T\{(T + k &- iEk A /fi + o^e'^'/*) |0) 

+ / <0| (cr,- + cr- 2 ) 2 e' w "' <r^ iHxrt {a + ki + < ) 2 |0> + 0(y 6 )\ . (S 15a) 


When measuring cr :, post-selecting measurement outcomes 
with a fixed number of total spin-up excitations (n) allows the 
dynamics (and therefore energy gaps) in different subspaces to 
be studied independently. We define projection operators, Il„, 
which denotes the projector onto the subspace with excitation 
number n. To study the single excitation subspace we extract 
the observable 


(Ht)\ nt (r\ n, |«A(0) 


C 2 y 2 


|*i> 


(S16) 


+ (k 2 \<r z i \k 2 ) 

+2 Re (e /(£i > ~ Eh >> tl% {k x \ cr] \k 2 )) 


Here, the signal oscillates with a frequency given by the en¬ 
ergy difference of the eigenmodes. For the initial state where 
tan(0,-) = y(A kl +A^ 2 ), the oscillations have amplitude pro¬ 
portional to 2A k ' A k l- (for k\ + k 2 ). The average signal from 
various spins is therefore maximized by Fourier-transforming 
the weighted sum 


MAt) = 2sign(A^A^ 2 )<n 1 cr5.n 1 ). (S17a) 

i 


Spectroscopic signals for two quasiparticle states 

The dynamics in higher excitation subspaces becomes more 
complicated because the quasiparticle modes that diagonalize 
Hxy in the single-excitation subspace are subject to the hard¬ 
core constraint cr] cr = 0. In a dilute system, i.e., when there 
are few excitations, J^-ctcrT/lV «: 1, the excitations become 
approximately independent. We can then construct a pertur¬ 


bation theory around the non-interacting eigenmodes without 
any hard-core constraint. 


Mapping spin operators to bosonic operators 

A formal way to construct a perturbation theory around 
the non-interacting eigenmodes is by applying the Holstein- 
Primakoff transformation [5] to obtain bosonic creation and 
annihilation operators b] and bj. 


(X; 


cr j 



-5 + b'b 


(SI 8a) 

(SI 8b) 

(S18c) 


where S = \ is the length of the spin (we keep the factor 2 S 
as it provides a convenient way to distinguish different pertur¬ 
bative contributions). The Holstein-Primakoff transformation 
conserves the bosonic commutation relations between spins 
on different sites [crj,cr^] = 0, i + j , and the onsite hard-core 
constraint crjcrj = 0. In linear spin-wave theory (LSWT) 
[6], one expands the spin Hamiltonian to leading order in the 
bosonic occupations, which amounts to neglecting the con¬ 


straint by setting 



1. The result is the approxi¬ 


mation 


H X y ~ Hlswt - 25 ^ J(b]bj + /z.c.) , (S19) 

<<i 



It is straightforward to diagonalize the linear spin-wave 
Hamiltonian //lswt to find its eigenenergies E k and eigen- 
modes fej = HjAjPj (the spin waves), which reproduce the 
exact physics in the single-excitation subspace. If there was 
no hard-core constraint, combinations of spin waves would 
also define the eigenstates in higher occupation subspaces. 
For low spin-wave densities, one can expect this constraint 
to only play a perturbative role, motivating the use of Hamil¬ 
tonian //lswt as starting point for a perturbative expansion in 
the spin-wave interactions. 

Perturbative treatment of spin-wave interactions 

One can obtain leading corrections due to hard-core inter¬ 
actions by expanding the square roots to next order in the spin- 
wave density ( b\bj) /2S, yielding 

H X y ~ //lswt + V , (S20) 

with V - -jT li jJij(PjP j bjbi + Pjb]bibiy The eigen¬ 
states of the unperturbed system are given by \ki,ko)sw — 
b\b\ 2 10) / sj 1 + <5a-[,/ t 2 , where |0) = is the spin-wave vac¬ 
uum and the 6 kuk2 takes care of correct normalisation. The 
perturbative shifts of the corresponding energies are in first 
order 

2(E kl +E h )M^ 

V klM = (k u k 2 \ V|*i,fc)sw =-p—-— • (S21) 

1 + Okik 2 

The strength of the interactions is determined by the geomet¬ 
rical overlap between different spin-wave modes, 

M h,k* = y A E A t 2A h A h ( S22 ) 

j 

and can be visualised as a contact interaction vertex between 
two incoming and two outgoing waves. In this picture, the 
factor 1/2 in the case k\ = k 2 can be understood from the 
fact that for indistinguishable waves there is only one outgoing 
channel, whereas there are two for distinguishable waves. 

Expected signal 

We can use these approximate energies to evaluate the ex¬ 
pected signal in the two-excitation manifold. To this end, we 
approximate the time evolution as 

e -/M \ kukl) x e-^wr+kvy^^sw + \6ift(k u k 2 ))). 

(S23) 

If we are interested only in the strongest frequencies, we 
can neglect corrections \6i//(k\ ,£2)) to the states \ki,k 2 ) = 
o+ a+ |0), which correspond to scattering into other spin- 
wave modes, and approximate 

& -iH X yt/ti | kukl ) x e -KE h +E k 2 + V h x 1 )tlh ^ ^ _ ( S24 ) 


We now define Ajt 2 = E kl - E kl . For any observable O 
that is diagonal in the tr Pauli matrices, the observed signal 
in the two-excitation subspace is 

m)\n 2 o z n 2 \m) * 

C 2 / [(k l ki\0 : \kih) + ( k 2 k 2 \O z | k 2 k 2 ) + 4{hk 2 \O z \hk 2 ) 

r* I — A E k | k-> F/7 1 .t' Vk 2 ,k 2 \ , . I .. . . 

+2 cos -—-_!_!- =LL t \ (kMO | k 2 k 2 ) 

a I A/(/.'i ,k 2 + VkiM ~ VkiJc 2 \ /j 7 1 .7 7 . 

+4cos I-—--AT- (k 1 k\ | O' \kik 2 ) 

A / A E k | k 2 F/„ l ,k 2 F^2 Jc 2 \ ,7 7. .7 7 . 

+4 cos -—-_LJ- =i± t 1 fk x k 2 \o I k 2 k 2 ) 

+ corrections from other states ] . (S25) 

This formula neglects all higher-order contributions in F, 
including scatterings into other than the initially prepared 
modes k\,k 2 . As long as perturbation theory is a good descrip¬ 
tion of the two-excitation subspace, we thus expect predomi¬ 
nantly three frequencies in the time-evolution of the magneti¬ 
zation signal, which lie at hv a = |A E klJc2 + V ki , kl - V ku k 2 1, hv b = 
AT/,., _k 2 + Vk 2 ,k 2 ~ F^ 2 ,fc 2 l» and hv c — |2A E kltk2 + F^, k \ ~ Ft 2j ^ 2 |. 
Indeed, this three-peak structure with v a + v b = v c is what we 
observe in the data for [k\ = l, k 2 = 7), as presented in Fig. 4 
of the main text. 

To derive simple quantitative estimates, we use sine- 
wave amplitudes for the eigenfunctions, A k = Jppp sin(k • 
jjpp]), yielding the interaction-shifted frequencies v 0 = (l — 
0.113)A E klJci /h, v b = (1 - 0.262)A E kiJc Jh, and v c = (2 - 
0375)AE klt k 2 /h for {k\ = \,k 2 =l). 

Summation of spin-spin correlation measurements 

In the two-particle subspace, we study weighted sums of 
two-spin projectors ’PE = | ( o~ + 1 j (tr 7 . + l), i.e., the prob¬ 
ability that both ions i and j are in the state It). The corre¬ 
sponding expectation values are 

^ {k x k 2 \ (crjy + 07 + cr. + l) | k 2 k A ) = 

(A k 'A k r + A k2 A k ') (A* 3 A* 4 + Ajviy) . (S26) 

For the component with oscillation frequency v c , we have 
\ < kiki\(cno~ + 07 + 0-5. + l)| k 2 k 2 ) = 4A k 'A k, A k2 A k2 , so that 
the signal is maximised by 

M 2a (t) = Yj si gn (A k ;A k 'A k2 A k2 ) (PE) . (S27) 

l<m 

For {k\,k 2 } = {1,7}, this choice simultaneously suppresses 
the expectation values associated to the other two oscillation 
frequencies, so that it allows one to cleanly extract a single 
one of the three expected Fourier components. 

Similarly, the component oscillating with v a has an ampli¬ 
tude 21A* 1 1 2 A k 'A kl + 2| A k > 1 2 A kl A k . 2 . If we use 

M 2h (t ) - J] ( si g n K'O + sign( a*A*'J) (P 1 ..) , (S28) 
i<j 
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Frequency (Hz) 

Figure SI: Spectroscopy of quasiparticle interactions. (a,b) 

Fourier spectra of the two-quasiparticle state | 4 ), obtained from 
summed signals (red curves), with M 2a {f) to enhance v c = v a + Vj 
(panel a) and M 2a (t) to increase signals at v ab (panel b). The mea¬ 
sured peak positions compare well with those of a simulated evolu¬ 
tion with longer times, yielding a narrower bandwidth (grey curves). 
However due to strong quasiparticle scattering a large number of 
frequencies components becomes apparent. The weakly interacting 
single-particle eigenmodes no longer provide a good basis for de¬ 
scribing the system. 


the signals due to this component and the one at v/, are maxi¬ 
mized while the one at v c gets suppressed. 

Validity of first-order perturbation theory: additional experimental 

data 

The quality of the first-order perturbation theory for given 
eigenstates can be estimated by considering how well a state 
\k\kf) remains localized in the two-excitation state space when 
taking the hard-core constraint into account. This amounts to 
acting on the vacuum |0) with b' k and computing the over¬ 


lap of the resulting state to the exact eigenstates. The state 
\k\ = 1,&2 = 7), constructed from single-particle states at the 
edge of the spectrum, for example, remains well localized in 
few exact eigenstates, and therefore qualitative predictions of 
simple first-order perturbation theory in the spin-wave interac¬ 
tion work quite well. Other combinations, on the other hand, 
such as |ki = 1,^2 = 4) (Fig. SI) show many other spectral 
features deviating from the three-peak structure predicted by 
the perturbation theory. This comes from a considerable over¬ 
lap with many of the true eigenstates in the two-excitation 
subspace, resulting in a diminished reliability of perturbation 
theory. 

Since the spin-wave scattering is proportional to the over¬ 
laps which decrease with reduced excitation density, 

estimates from perturbation theory will improve with increas¬ 
ing system size (if the excitation number is kept fixed). For 
the small system considered, the overlaps are quite large, ly¬ 
ing in the range A» 0.1 - 0.25, which explains the lack 
of quantitative reliability of non-linear spin-wave theory. For 
N -r 00 , however, corrections such as go to 0 as 1 /N, 

and the states |kik2)sw P rov id e accurate approximations to the 
true eigenstates of the two-excitation subspace of Hxy- 
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